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ABSTRACT 

A selection of statistically stable (robust) algorithms for data variance calculating has been 
made. Their properties have been analyzed via computer simulation. These algorithms would 
be useful if adopted in radio astronomy observations in the presence of strong sporadic radio fre- 
quency interference (RFI). Several observational results have been presented here to demonstrate 
the effectiveness of these algorithms in RFI mitigation. 

Subject headings: methods: miscellaneous — methods: statistical 



1. Introduction 

Cosmical radio emissions received by radio tele- 
scopes are noise-like signals which are character- 
ized by a normal (Gaussian) probability distribu- 
tion function A/'(0, (7^^^), i.e., with zero mean and 
variance cr?- • Background radio emission and ra- 
dio receivers also produce normal noise A/'(0, cFsys)- 
Given a large number of n samples Xi^i = 0..n 
of normally distributed data with zero mean, a 
classical statistical procedure for the estimation of 



variance cr" is 



n 



(1) 



The time-honoured estimate (1) may be ob- 
tained by usin^ the Maximum Likelihood method 
( van der Waerden1ll969. ) and is a minimum vari- 
ance unbiased estimate. Traditional radiometrical 
techniques implement this procedure in analogue 
or digital form in conventional radiometers: the 

^ is obtained after square- 
ass filtering (averaging) 
. This statistical infer- 



estimate of a?- + 
law detection and 



^ sys 

low- 



( Rohlfs fc Wilson M2QQ, 



ence is based upon an essential assumption about 
the probability distribution of the data "signal 
noise + systen noise": the distribution must be 
normal. As in many cases in applied statistics 
this assumption is not always valid for radio as- 



tronomy. Most of the data may correspond to 
the model with an assumed normal distribution, 
but there are also a number of outliers^ atypi- 
cal data which stand out from the bulk of the 
data. As a result there is an "approximately" 
normal distribution that gives rise to outliers, i.e., 
the distribution has a normal shape in the cen- 
tral region but has tails that are heavier than 
those of a normal distribution. If such approxi- 
mate normality were to hold, the results of using 
a normal distribution theory will nevertheless not 
hold approximately. In the presence of heavy tails 
the estimate (1) based on the maximum likelihood 
principle is no longer the best and may have unac- 
ceptably low statistical efficiency (large variance) 
and very large bias. These phenomena in applied 
statistics have been well known from the time of 
Gauss, Newcomb and Eddington. 

Radio frequency interference (RFI) creates a 
situation in radio a stronomy where such out- 
liers arise ( Fridman fc BaanJ I2QQII : [k 



s^esteven 



20071 ). Industry, ground and satellite commu- 
nications, ground and airborne radar, power 
lines, radio and TV stations, etc. produce all 
kinds of additional noise which penetrate into 
extremely sensitive radio astronomy receivers. 
Radio observatory computers are themselves of- 
ten the sources of RFI. Many authors demon- 
strated the effectiveness of different RFI mit- 
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igation methods (lEllingson fc Hampson 2QQ3I: 



Baan, Fridman & Millenaar ' '2004'; "Fisher et al. 



2005; Mitcheh, Robertson & Sault 2005; Poul sen. Jeffs 
20051 : IZhang et al. 1 120051 : iKesteven I i2007) in ob- 
servations with single dishes and radio interferom- 
eters. Adaptive noise cancellatio n (ANC) with ref- 

erenc e antenna was described in (jBarnbaum fc Bradlev 
19981 ). Post detector and post correlation pro- 
cess ing: using; reference a nte nnas was suggested 
in (|Briggs.Bell fc Kesteven I I2OO0I ). Higher or- 
der statistics ( HQS) analysis was proposed in 
( Fridman 200l[ ). Principal results were achieved 
with blanking and reference antennas. The high 
expectations of new projects (LOFAR, ATA, 
SKA) are based on the implementation of spatial 
processing - nulling in the direction of RFI. But 
very often RFI do not obey the model assump- 
tions which guarantee the success of this kind of 
processing. They are variable, sporadic and are 
not fully coherent at the s i tes of multi-element 
syste ms ([Thompson I 120031 : I Jeffs. Li fc Warnick " 
20051) . 



In this paper a new view on blanking is pro- 
posed based on statistical methods granting stable 
estimates of variance in the presence of outliers. 
These methods would be useful in radio astron- 
omy RFI mitigation, including real-time process- 
ing. The aim of this paper is to compile a selection 
of methods of stable variance estimates suitable 
for RFI mitigation and to study their properties 
via computer simulation. Some examples obtained 
during real observations are also given. 

The typical appearances of RFI is of strong 
impulse-like bursts in both temporal and spectral 
domains. They behave as outliers in the data with 
normal distribution and the data received may be 
characterized by the contaminated normal distri- 
bution 

F{x) = (1 - e)N'iO, asys) + eFRFi{x), (2) 

where A/'(0, dsys) is the "clean" probability dis- 
tribution, FiiFi{x) is the usually unknown dis- 
tribution of RFI, e characterizes the fraction of 
Frfi{x) in the total < e < 1. The theory 

of ro b ustness p ioneered by (jTukev Ill960l : Huber 



1963 iHampel I Il971) gives recommendations on 
obtaining statistically stable estimates of parame- 
ters in the situation similar to (2). The estimate is 
called statistically stable (robust) if slight changes 
in distribution have a relatively small effect on its 



value, or, in other words, robustness means insen- 
_sitivity of a statitical procedure to small deviations 
^J&mk^umptions about the normality of data's 

distribution. 

Several approaches to developing robust esti- 
mates exist: M-estim ates based on a generalized 
maximum likelihood (|Huber Ijgool. L -estimates 



based on order statistics (iDavid fc Nagaraia 



2003_) , i?-estimates based on ranks (jHettmansbperger 
1984f ) an d estimates ba sed on nonparametrical 
statistics (|Sheskin Il2000[ ). 

2. Stable estimates of variance 

A selection of stable (robust) estimates of vari- 
ance will be described in this section. There are 
different ways of judging the stability of an es- 
timate. Given a distribution F{x) and a popu- 
lation X representing this distribution, an esti- 
mate T(X, F) is the functional depending on F{x). 
If this estimate is relatively unaffected by small 
changes in F, i. e., the functional is continu- 
ous with respect to F, then the estimate is said 
to have qualitative robustness. For a distribution 
F{x) contaminated by a distribution G{x) with 
probability e 

F,{x) = (1 - e)F{x) + eG{x) (3) 

the influence function (jHampel et al. Ill986h /F(T, F) 
is widely used in theoretical analysis 



IF{x) lim 



T(F,,e) - T{F,) 



(4) 



The /F(T, F) is the relative influence of x on T{F) 
having the value x with probability e. The esti- 
mate T{F) is said to have infinitesimal robustness 
if IF{x) is bounded. The IF of classical both 
arithmetical mean and variance are not bounded. 
A simple flnite-sample version of IF called em- 
pirical influence function (or sensitivity function) 
exists. It is constructed in the following way. 
Suppose there is an estimate T^-i of a sample 
(xi, Xn-i). A new n — th sample x is added 
and Tn(xi, ...,Xn-i,x) is considered as a function 
of X. For the estimate (1) the sensitivity function 
is 

IFi.) = ^ (5) 

Another useful characteristic of robustness is the 
breakdown point: the minimum value of e for which 
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an estimate goes to infinity as x grows. The break- 
down point of conventional mean and variance is 
0. The breakdown point of median is 0.5, i.e., up 
to 50% of data may be outhers and they do not 
alter the correct value of the median. 

Taking into consideration our practical goal to 
use estimates of variance in observations, two mea- 
sures of stability will be exploited for the following 
statistical procedures. 

1. The relative empirical influence function 
with the following modification: for a sample 
length N the number of additional outliers with 
an amplitude Xrfj will be M = [eiV], where [q] 
means the greatest integer part of q. The estimate 
calculated for a contaminated sample is compared 
with the estimate correspondimg to a "clean" sam- 
ple: 

REIF{ex) = TN+M{F.,e)-TN+M{F,,e=0) ^ 

The relative empirical influence function will be 
calculated during computer simulation for the es- 
timates described below. 

2. Stability against outliers is achieved at the 
expense of the effectiveness of an estimate. In the 
absence of outliers the standard deviation of a ro- 
bustly estimated variance is, as a rule, larger than 
that of a simple estimate (1). To characterize this 
loss the parameter LOSS will be used: 



SNRt^ 



SNRO 



LOSS = SNRtj,/SNRO, 

rms[TN{cr)] 
TOiv(oM^Acr) - TO^) 
rms[TON{o-)] 



(7) 



For a small increase of Act « a two signal- 
to- noise ratios SNR are compared: SNRt^ of a 
given estimate of variance TN{a) and SNRO of the 
estimate (1) with the best potential effectiveness. 

The relative empirical influence function and 
loss are calculated via computer simulations (cr = 
l,n = 10^) for the estimates presented below. 

2.1. Variance of the trimmed data 

Let be a random sample and let 

^(1) ^ ^(2) ^ ••• ^ X(n) be the observations sorted 



in ascending order. The ith largest value is 
called the ith order statistic. Let 7 denote the 
chosen amount of trimming, < 7 < 0.5 and 
k = [777,]. The sample trimmed variance is com- 
puted by removing the k largest and k smallest 
data and using the values that remain: 



N-k 

1 — AT- r»7 / ^ y^n l^trim) 
n=k 



N -2k 



(8) 



N-k 



-9k ^ 



N-2k 



n=k 



where jitrim is the sample mean of the trimmed 
data. Trimming lessens the variance of data and 
the coefficient Ktrim makes Ti the consistent esti- 
mator for data with normal distribution. Table 1 
gives the values of Ktrim for different 7. 

Table 1: Consistency factor and LOSS as func- 
tions of 7 for trimming 



7 


0.005 


0.01 


0.025 


0.05 


0.1 


Ktrim 


L085 


1.147 


1.32 


1.6 


2.28 


LOSS 


0.98 


0.95 


0.90 


0.86 


0.83 



The third row in this table gives the values of 
LOSS for different 7. 

Fig. 1 shows the relative empirical influence func- 
tions for untrimmed data (Fig. la) and for the 
variances computed for trimmed data with 7 = 
0.05 (Fig. lb). The parameter e ("eps") indicates 
the percentage of outliers in the total volume of 
data. 

2.2. Winsorized sample variance 

A sample xi,....,^^ is sorted in ascending or- 
der. For the chosen < 7 < 0.5 and k = [yn] 
winsorization of the sorted data consists of setting 

{^(/c+l), ^/ ^(i) < ^(/c+1) 
if ^(/c+l) < ^(i) < ^(n-/c) 
X{n-k), ^(i) > ^(n-/c) 

The winsorized sample mean is /i^ = ^ SiLi 
and the winsorized sample variance is 



To 



1=1 



(10) 



Table 2 gives values of the consistency factor 
Kujinsor for different 7. 
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Table 2: Consistency factor and LOSS as func- 
tions of 7 for winsorization 



7 


0.005 


0.01 


0.025 


0.05 


0.1 


Kwinsor 


1.019 


1.041 


1.103 


1.227 


1.54 


LOSS 


0.99 


0.98 


0.96 


0.88 


0.82 



The third row of this table gives the values of 
LOSS for different 7. 

Fig. Ic shows the relative empirical influence func- 
tions for the variances computed for winsorized 
data with 7 = 0.1. The parameter e ("eps") indi- 
cates the percentage of outliers in the total volume 
of data. 



2.3. 



Median absolute deviation 



This estimate for sorted data < X(2) < ••• < 
X(^) is defined by 

T3 = 1.483 X medi<i<ri{\xi — med{xi)\}^ (11) 

where 

med = 0.5(x(^) + n = 2m, 

med = n = 2m + 1 

The breakdown point for this estimate is 0.5, i.e., 
almost half the data may be contaminated by out- 
liers. But the effectiveness of T3 is much lower 
than for the Ti and T2: the LOSS = 0.6. 



2.4. Interquartile range 

This estimate is defined by 

T4 = [(xl_,-x,)/1.35]^ 



(12) 



where q = 0.25 and for distribution V Xq satisfies 
V{x < Xq) = q. The breakdown point for T4 is 
0.25 and the LOSS = 0.6. Fig. li shows the rela- 
tive empirical influence functions for the variances 
computed with T4. 

2.5. Median of pairwise averaged squares 

This estimate is motivated by the Hodges- 
Lehmann estimate of location 



, X 2 I X o , 

med[ 1 < ^ < J < 



which for variance is 

T5 = 1.483 X medi<i<jn[{x^, + x^^)/2] (13) 



Fig. Id shows the relative empirical influence 
functions for the variances computed with T5 . 
The parameter e ("eps") indicates the percentage 
of outliers in the total volume of data. The LOSS 
for T5 is 0.78. 

2.6. Qn estimate 



Th is estimate is proposed in (jRousseuw fc Croux 
1993h and requires fewer operations than T5, but 
is quite effective. It combines the ideas of the 
Hodges-Lehmann e stimate and the Gini estimate 
(jKendal fc Stuart 1^967^ : 



Te = 2.2219{\xi-Xj\,i<j} 



(fc), 



(14) 



where k = (^) and h = [n/2] + 1, A: ^ (^) /4. This 
estimate is the kth order statistic of the (2) in- 
terpoint distances. Fig. le shows the relative em- 
pirical influence functions for variances computed 
with Te. The LOSS = 0.88. 

2.7. Bi weight variance 

First the following auxiliary values are calcu- 
lated 

Xi - M 

J- i 



9 X MAD' 

where M and MAD are the median mean and 
median absolute deviation, respectively. Then the 
coefficients are defined as 



aj 



_ r 1, if\Y,\<i 



0, if\Yi\>l 
in which case the estimate is 



Fig. If shows the relative empirical influence 
functions for variances computed with T7. The 
LOSS = 0.92. 



2.8. Bend Mid variance 

This estini ate of variance is described in 
(| Wilcox 112001 ). Set f3 = 0.1 and m = [(1 - f3)n + 
0.5]. Let Wi = \xi — M\ = 1, ...n, and are 
these numbers sorted in ascending order, M is the 
median mean of data x. Then the estimate of 1 — /3 
quantile of the distribution of W(^i^ is cJ^ = • 
Now set 

_ Xj-M 

^ i — 
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and 



1, if\Y^\<l 
0, if\Yi\>l 



in which case the bend midvariance is 



where 



= max[— 1, min(l, t)]. 



(16) 



Fig. Ig shows the relative empirical influence 
functions for variances computed with Tg. The 
LOSS = 0.96. 

2.9. Estimate with exponential weighting 

This estimate is described in ( Shurvgin 



Mean is estimated as a solution of the equation 

n 

Jl;: ^(xi-/i;)e-^^/'^ = (17) 

i=l 

and variance is defined as a solution of the equa- 
tion 



-Qi/i 



0, 



(18) 

where qi = {xi — /iP)^/cr^. Fig. Ih shows the 
relative empirical influence functions for variances 
computed with Tg. The LOSS = 0.98. 

D. A. Lax performed a Monte Carlo study of 
more than 150 variance estimators. Seventeen 
of these estimators were sele cted as bei ng either 
promising or commonly used (|Lax Ill985[ ). The re- 
sults of our co mputer simu lation are in agreement 
with this test ( Lax especially concerning 

the high efficiency of estimate Ti and T7 and the 
low efficiency of T3. 

The choice of a particular estimate depends on 
the type and intensity of RFI, the type of obser- 
vations and the method of implementation (hard- 
ware or software). and Tg are the best esti- 
mates from the point of view of LOSS. T3 and 
T5 remove outliers in a most effective way (high 
value of the breakdown point). Number sorting 
and permutations of pairwise measurements which 
are necessary in several algorithms require more 
computational time and computer memory. 



3. Computer simulations 

This section presents computer simulations of 
RFI mitigation using algorithms which were em- 
ployed in observations (Section 4). Two esti- 
mators from the the previous section have been 
chosen: winsorization and exponential weighting. 
The reasons for this choice are as follows. 

All algorithms described in Section 2 give esti- 
mates of variance, i.e., they work as total power 
detectors (TPD) in radiotechnical terminology. 
They can be applied in single dish observations 
both in continuum and in spectral observations. 
In the same manner they can substitute for TPDs 
which are already installed in existing radio tele- 
scopes. Nowadays it is practically impossible for 
technical and organizational reasons. In future ra- 
dio relescopes may be equipped with some RFI 
mitigation techniques and these variance estima- 
tors could then be implemented in hardware or 
software shape. But at the present time any exper- 
iment with RFI mitigation at existing ratio tele- 
scopes must take the technical constraints of im- 
plementation into account. Therefore only those 
estimators from Section 2 were chosen which could 
provide not only estimates of variance but also 
"clean" data which could be applied further to 
TPD or to the correlator already in use in a radio 
telescope backend. 

3.1. Winsorization followed by total power 
detectors 

Winsorization with the parameter 7 = 0.05 
was chosen mainly because impulse-like RFI in the 
temporal domain was predominant during obser- 
vations made at Effelsberg presented in Fig. 8. It 
was strong and sparse. The percentage of RFI in 
the whole volume of data was not large (< 5%). 
Both winsorization and trimming can provide raw 
data without outliers but trimming inconveniently 
reduces the number of samples because of this im- 
plementation problem: "timing" is distorted, i. e., 
there are gaps (not zeros) in the presumed un- 
broken flow of data and the radio telescope back- 
end has no information about the location of these 
gaps. The problem is easily solved in software dur- 
ing simulations but may be critical in real obser- 
vations. Winsorization unlike trimming preserves 
the total number of samples. 

The block diagram of computer simulation pre- 



5 



sented in Fig. 2 shows the structure of the al- 
gorithm using winsorization. The input signal is 
the mixture of three signals: system noise with 
normal distribution (a = 0.5), signal noise, also 
with normal distribution (cr = 0.05) and impulse- 
like interference imitating RFI. The signal noise 
is switched "on" and "off" emulating "on-source" 
or "off-source" position of the antenna main beam, 
respectively. Poisson distribution (A = 0.004) gov- 
erns the appearance of RFI and the lognormal 
distribution (mean=12, standard deviation=6) de- 
termines the random amplitudes of impulses. In- 
put signal samples are stored in the buffer and 
sorted to create order statistics. The length of 
these statistics is equal to N = 100. The win- 
sorization block provides both variance estimates 
for each of the samples and also the winsorized 
samples which are applied to the external TPD (as 
in observations presented in subsection 4. 1). The 
loop with M = 100 cycles corresponds to post- 
detector averaging: the estimates of variance after 
the winsorization block are accumulated during M 
cycles and the mean calculated for these M cycles 
represents the simulation result. 

Fig. 3 shows the results of computer simula- 
tions : 

a) input noise with normal distribution, /i = 
0.0, cr = 0.5, no interference; b) total power detec- 
tor output, each point in this figure corresponds to 
squaring and averaging of TV x M = 10^ samples 
in figure a); there are two steps, "up" at point 
#100 and "down" at point #200 corresponding 
to the increase of a from value 0.5 to the value 
a + Act, Act = 0.05 ( "on-source" and "off-source" 
position of antenna beam); 

c) interferences in the form of random impulses 
are added to the noise a); 

d) total power detector output with the input sig- 
nal c); 

e) total power detector output with the input sig- 
nal c) and preliminary winsorization. The differ- 
ence between the scales in the vertical axes in d) 
and e) demonstrates the effect of interference mit- 
igation. 

3.2. Exponential weighting 

Exponential weighting also provides not only 
estimates of variance but data with suppressed 
outliers. This estimator is also of superior effec- 
tiveness. And unlike winsorization, it does not 



sort input data. The influence of sorting on phase 
information can be detrimental to radiointerfer- 
ometic and pulsar (de-dispersion procedure) ob- 
servations. Impulse-like RFI in the frequency do- 
main was prevalent during the pulsar and image 
synthesis observations at WSRT (subsection 4. 2 
and 4. 3). It is to demonstrate exponential weight- 
ing in the frequency domain which is the aim of 
this computer simulation. 

The block diagram presented in Fig. 4 shows 
the structure of the algorithm using exponential 
weighting in the spectral domain. The input sig- 
nal is a mixture of four signals: system noise with 
normal distribution (a = 1.0) filtered to emulate 
the WSRT backend low-pass filters; two noise-like 
signals with normal distribution imitating emis- 
sion and absorption spectral lines, respectively, 
and interference waveforms. A randomly phase- 
modulated sinusoidal carrier is used as the RFI 
signal, the index of modulation is equal to 7r/2 
(binary phase modulation). Poisson distribution 
(A = 0.3) determines the moments of phase jumps. 

There are two loops in the algorithm: the inner 
loop (m = 1, ..M) which is responisble for the stor- 
age of M instantaneous spectra of the input signal 
for RFI mitigation processing, and the outer loop 
(/ = 1, ..L) which repeats consecutively L identical 
stages of processing inside the inner loop, showing 
the dynamics of the running spectrum. Interfer- 
ence is intermittently switched "on" and "off" in- 
side the inner loop, the percentage in the total 
averaging interval is less than 7%. 

N samples of the input signal are Fourier- 
transformed with forward FFT, N = 512. The 
inner loop with the counter (m = 1, ..M) provides 
M complex spectra (each having a length equal to 
N) which are stored in the buffer. So there are 
two sets (real and imaginary) of data in the buffer 
numbered as m = 1, ..M for each of the = 512 
spectral channels. L power spectra of the input 
signal averaged on M instantaneous spectra are 
calculated in the outer loop (/ = 1,..L). Fig. 5a 
shows the three-dimensional time-frequency pre- 
sentation of the sequence of L = 50 averaged 
"dirty" spectra. 

Equation (18) is solved separately for each of 
2 X N real and imaginary sets of data providing 

estimates of aj^ 

,reaii'^) ^\i,imag^^) ~ variances 
in i = 1, ..AT spectral channels calculated using M 
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samples of real and imaginary components of the 
instantaneous spectra. Random values in the real 
and imaginary part of the instantaneous spectrum 
are independent, so the sum of these estimates 
in each spectral channel gives the total estimated 
power spectrum for each spectral channel. This 
"clean" power spectrum is represented in Fig. 5b 
for L = 50 time intervals. In real observations 
each oil — th time interval is equal to N x M x /Si 
seconds, where At = 1/2 A/ is the input signal 
sample interval. A/ is the bandwidth of the input 
signal. 

In appHcations it is practical to find the solution to 
(18) when fir = using the approach of stochastic 
approximation: 



^ m 



^jexpi- 



(19) 

m = 1, ..M. Fig. 5c and 5d show the re- 
sult of averaging L spectra in Fig. 5a and 5b 
for the "dirty" and "clean" power spectra, respec- 
tively. The frequencies of interferences do not co- 
incide with the frequencies of spectral lines . It can 
therefore be seen that the components with nor- 
mal distribution (system noise and spectral lines) 
are untouched by this RFI mitigation procedure. 
The sequence of pictures in Fig. 6 similar to Fig. 
5 demonstrate the result of computer simulation 
when the frequencies of interferences coincide with 
the frequencies of spectral lines . In this case we 
again see that the restoration of the spectra is sat- 
isfactory. 

Figures 5 and 6 demonstrate how the algorithm 
works with spectra. However experimental con- 
straints did not require the provision of estimates 
of variance or power spectrum but did required the 
"clean" signal in the temporal domain similar to 
the input signal (the level of the signal, the band- 
width). Therefore running estimates of "clean" 
and "dirty" power spectra are used for exponential 
weighting of the running complex instantaneous 
spectra: each complex value of an instantaneous 
spectrum being multiplied by exj9(— P(z)/3<j^(i)), 
where P{i\i = 1, ..A/" is the input spectral vari- 
ance (power spectrum) in the channel i, g^(i) is 
the estimate of the "quiescent" power spectrum in 
the channel i found from Eq. (18). A correspond- 
ing delay must be introduced because of the time 
required to calculate all g^(i). Then, after the 



backward FFT, the signals in the temporal do- 
main can be applied to TPDs or correlators. This 
auxiliary output was used to supply the "clean" 
raw data in the observations described in Sections 
4. 2 and 4. 3. 

Figure 7 shows the impact of RFI mitigation 
on cross-correlation (radiointerferometric observa- 
tions). The algorithm in the block diagram shown 
in Fig. 4 was applied to two signals with the ad- 
ditional coherent Gaussian component (gs = 0.1) 
emulating the noise from a radio source received 
at both sites of the radio interferometer. Sporadic 
interference with the duty cycle equal to 0.2 was 
generated as frequency-modulated carrier and was 
identical at both sites which is the worst case sce- 
nario: 100% correlated RFI. In reality the impact 
of RFI in radio interferometers is considerably re- 
duced due to fringe stopping and delay tracking 



procedures (^Thompson 1 1198 



Figures 7a and 7b demonstrate three-dimensional 
time-frequency presentations of the input "dirty" 
and "clean" power spectra on one site, respec- 
tively. Deep "troughs" can appear in the spec- 
trum in Fig. 7b at the places of interferences due 
to the exponential effect of weighting: when in- 
terference is strong the algorithm works similarly 
to the "thresholding and blanking" algorithm but 
more smoothly and without the a priori knowl- 
edge necessary for the positioning of the threshold 
level. Figures 7c and 7d give the averaged nor- 
malized cross-correlation functions corresponding 
to Fig. 7a and 7b, without and with RFI mitiga- 
tion, respectively. The central parts of the cross- 
correlation functions between channels #400 and 
#600 are shown. The scale of the vertical axe in 
Fig. 7c is significantly larger than that in Fig. 
7d thus showing strong excessive cross-correlation 
due to RFI. 

It is necessary to add two further comments 
here for the purposes of discussion. 

1. The estimators of variance described in Sec- 
tion 2 and those tested in Section 3 are essentially 
nonlinear procedures. For example, one can, at 
least theoretically, imagine a huge outlier which 
renders the conventional estimator (1) completely 
non-functional whereas these estimators work per- 
fectly. Radioastronomers often ask what is the 
value of RFI suppression (usually in dB) when an 
RFI mitigation algorithm is applied? Example in 
Fig. 3 demonstrates that it is not always possi- 
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ble to answer this question correctly. The larger 
the outlier's amplitudes the more effectively they 
are deleted. If the amplitudes of impulses in Fig. 
3c are 10 times larger the result in Fig. 3e will 
be the same or even slightly improved. This is 
also valid for the "thresholding and blanking" al- 
gorithm which has been succesfully used in several 
articles on RFI mitigation. That is why it is bet- 
ter to judge an RFI mitigation procedure in each 
paticular case in a combined manner, that is look- 
ing at the "signal-of-interest" distortions by both 
the RFI residuals and the procedure itself. 

2. The estimators of Section 2 are intended for 
applications when sporadic, impulse- like interfer- 
ence disturbs observations. However interference 
is often continuous, persistent and practically con- 
stant during the averaging interval. For example, 
for M = 100,7V = 512, A/ = 20MHz this inter- 
val is equal to 1.28 x 10~^ sec which is rather a 
short interval for several RFI to change their am- 
plitudes. The estimators of Section 2 will not yield 
any benefits in this case. But the same fact that 
RFI is quasi-constant at a reasonably chosen aver- 
aging interval can help to decouple the Gaussian 
component with a "quiescent" value of a and RFI 
with ease. A simple procedure can do this. The 
power spectrum after FFT is usually calculated as 

M 

P{f) = l/M Y.lMsifU? + Im[s{f)m?} 

m=l 

(20) 

The sums of RFI and "useful" Gaussian noise 
in each spectral channel are squared and aver- 
aged. Strong RFI dominates in this case over the 
noise. Another algorithm is proposed here: to sep- 
arately estimate the variances of real and imagi- 
nary parts of the spectrum and then add them to 
obtain the total variance in each spectral chan- 
nel, i.e., the power spectrum sought. The vari- 
ances of real and imaginary parts of the complex 
spectrum must be calculated taking into consider- 
ation that the mean value now is no longer equal 
to zero because of the presence of RFI. So the 
general formula for sample variance must be ap- 
plied: for a random value x, the sample variance 
is var{x) =cr^ = '^/MY^^^-^{xm — x)'^ which for 
the power spectrum at frequency / is 



P{f) = var[Re{s{f)] + var[Im{s{f)] = (21) 

M 

l/M ^ Re[s{f)]m - R4^f)]V + 

m=l 
M 

l/M^/m[s(/)U-/MK/)]}', 

m=l 

M 

= l/M ^ Re[s{f)]m, 

m=l 
M 

I^f)] = l/M J2 ImHf)]m 

m=l 

In this case RFI is eliminated due to its con- 
stant value for all M samples of the instantaneous 
spectra s{f)]m' It may also be useful for the mit- 
igation of weak but persistent RFI whose detri- 
mental impact is revealed only after lengthy aver- 
aging. This is a conjecture which has been proved 
in computer simulations but which, of course, re- 
quires experimental confirmation. 

4. Examples from observations 

4.1. Observations in continuum at Effels- 
berg radio telescope 

Radio source 1448+762 was observed in con- 
tinuum at the Effelsberg radio telescope. The re- 
ceiver output was split on two channels. The sig- 
nals from one channel were processed by a FPGA 
processor (Altera Stratix S80) and then sent to 
the total power detector. The signals from the 
second channel were applied straight to the total 
power detector. The bandwidth of the signals ap- 
plied to the total power detector was equal to 20 
MHz. These channels provided two radio telescope 
outputs: one output with RFI mitigation and an- 
other without. Pairs of radio source scans from 
these outputs were made simultaneously. Because 
the total power detector was an integral part of 
the radio telescope backend equipment, the algo- 
rithm implemented in FPGA only processed the 
IF (intermediate frequency) signal with the aim 
of "cleaning" it of RFI. The analogue input sig- 
nal was digitized in 12-bit ADC with 40Msamples 
speed, then processed in FPGA, transformed back 
into analogue form and applied to the total power 
detector. Winsorization of the signal in temporal 
domain (7 = 0.05) was implemented. 
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Fig. 8 shows eight scans of the radio source, 
each scan represented by two panels: the top panel 
- scan with RFI mitigation, the lower panel - with- 
out RFI mitigation. Fig. 9 displays radio images 
of the source built using scans similar to those in 
Fig. 8: the left panel - without RFI mitigation, 
the right panel - with RFI mitigation. 

4.2. Pulsar observations at WSRT 

A new pulsar machine PUMA-2 has been in- 
stalled at the Westerbork synthesis radio telescope 
(WSRT). The radio telescope works in tied-array 
mode in which all 14 signals from antennas are 
added in phase , i.e., there is one output (in reality 
with two polarizations) as for a single dish. The 
20 MHz-baseband signals from each of the eight 
frequency channels of WSRT are digitized (8 bit) 
and stored in the mass storage system which has 
sufficient hard disk capacity to support at least 24 
hours of continuous observations. Signal process- 
ing can therefore be undertaken off-line. 

In our experiment a block of data recorded dur- 
ing 10 sec, 40 X 10^ eight-bit samples/sec, was 
used. All RFI mitigation processing and the to- 
tal power detector (TPD) were realized entirely 
in software during off-line processing. The esti- 
mate of variance with exponential weighting was 
used. Fig. 10 displays the results of processing. 
Upper row, left panel: TPD outputs for two po- 
larizations calculated from the raw data with RFI, 
right panel: TPD outputs, RFI removed. Middle 
row, left panel: example of a time fragment of the 
power spectrum with RFI; right panel : the same 
time fragment, RFI removed. Lower row shows 
pulsar profiles after de-dispersion and folding at 
the pulsar period, left panel: pulsar profile aver- 
aged over 10 sec from both polarizations of raw 
data, middle panel: similar pulsar proffie, RFI re- 
moved, right panel: pulsar profile restored with 
observational data obtained at 1420 MHz without 
any RFI which is put here for comparison with the 
profile in the middle panel. 

4.3. Radio image synthesis at WSRT 

Radio source DA240 was observed at WSRT at 
a central frequency of 357 MHz with the band- 
width equal to 20 MHz in the presence of strong 

RFI. The RFI mitigation s ystem (RFIMS) was 

used for real-time processing (iBaan. Fridman fc Millenaar 



20041 ). Analogue baseband signals were digi- 
tized (12bit ADC, 40 Msamples/sec), processed 
in FPGA (Altera StratixS80) and transformed 
back to analogue form for subsequent processing 
in the WSRT correlator. The algorithm used for 
the removal of RFI was similar to the algorithm 
with exponential weighting, except that the vari- 
ance was not used as an output, instead "cleaned" 
, exponentially weighted signals were applied to 
the correlator. 

Radio images of the source DA240 are shown 
in Fig. 11. Upper row, left panel: image with- 
out RFI mitigation; right panel: image with RFI 
mitigation. Lower row: central parts of the image 
presented in the same order. The stretched form 
of the synthesized images is explained by the fact 
that the observations lasted 8 hours, instead of a 
full 12 hour aperture synthesis cycle. 

RFI mitigation was implemented on each of 
14 radio telescopes at WSRT. This may give rise 
to some distortions because of the difference of 
equipment characteristics at the different anten- 
nas. Other observations were made specifically to 
judge the "toxicity" of the RFI mitigation proce- 
dure. Two pairs of radio images of radio source 
4C34.47 were synthesized: one was observed at 
1420 MHz with the RFI mitigation system and 
without it and another at 345 MHz, also with the 
RFI mitigation system and without it. The ffist 
pair of images served as a test for "toxicity" while 
the second pair demonstrated the effectiveness of 
RFI mitigation. Both pairwise observational data 
were stored during simultaneous observations. 

Fig. 12 displays the results of the image synthe- 
sis. Upper row, left panel: central frequency 345 
MHz, without the RFI mitigation system; right 
panel: central frequency 345 MHz, with the RFI 
mitigation system. Notice the difference of inten- 
sity levels in the figures. Lower row, left panel: 
central frequency 1420 MHz without the RFI miti- 
gation system, right panel: central frequency 1420 
MHz with the RFI mitigation system. After sub- 
tracting one image from the other the rms noise 
is less than 0.7 mJy/beam which signifies a good 
similarity in the images. 

5. Conclusions 

1. Statistical analysis of raw data with the 
finest available time and frequency resolution can 
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help during observations in an RFI contaminated cessing of 4C34.47 data, 

envinronment. Estimates of variance are an im- 
portant part of both classical and robust statistics. 
Growing concern about RFI pollution should per- 
suade the radio astronomy community to pay more 
attention to a variety of algorithms developed in 
the realm of robust statistics. Tradional radiome- 
ters with simple square-law detectors, spectrum 
analyzers and correlators must be equipped with 
these tools. This framework of robust estimates 
puts the successfully tested blanking of RFI on a 
more stable foundation. 

2. Statistically faithful, robust estimates of 
variance are especially appropriate for application 
in an impulse- like strong RFI environment. RFI 
is effectively suppressed and the accompanying 
loss in the signal-to- noise ratio is tolerable. Ta- 
ble 3 gives a summary of available RFI mitigation 
methods which have been proposed during the last 
years. Some of these have been tested in real ob- 
servations. Blanking showed truly good results in 
the case of impulse-like strong RFI. There are slots 
in the Table 3 where "removing or blanking"" is 
referred to. The aforementioned robust algorithms 
can be applied in these particular situations. 

3. The results produced from the observations 
at the Effelsberg radio telescope (image synthesis 
by beam scanning) and at WSRT (Earth rotation 
aperture synthesis) included in this article demon- 
strate the usefulness of this kind of RFI mitigation 
both in real-time and off-line (pulsars). 

4. The choice of a particular algorithm depends 
on the type and intensity of RFI. The proportion of 
RFI presence in data is also important. The type 
of implementation may determine the choice: off- 
line or real-time. All options are open now. Mass 
storage systems allow data to be processed off- 
line on powerful computers. Nowadays both exist- 
ing radio telesopes future projects (LOFAR, ATA, 
SKA) will generate such huge amounts of data 
that real-time processing is vital: DSP, FPGA or 
supercomputers are possible solutions. Therefore 
the creation of " robust" radiometers, spectrum 
analyzers and correlators is an urgent necessity. 

I am grateful to Juergen Neidhoefer and Ernst 
Fuerst for their help during observations and data 
processing at the Effelsberg radio telescope, to Ben 
Stappers and Ramesh Karuppusamy for providing 
the pulsar data and Subhashis Roy for AIPS pro- 
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Table 3: The application of RFI mitigation methods 



Type of radio 
telescope 


Type of RFI, 
intensity 


Type of RFI, 
structure 


Observations in continuum 


Spectral observations 


Single 
dish 


weak 
RFI 


impulse-like 


tolerable 


tolerable 


narrow-band 


Removing (blanking) 
in spectral domain 


ANC with a reference channel 
and post-TPD subtraction; 
HOS analysis 


wideband 


ANC with a reference channel 
and post-TPD estimation 
& subtraction 


ANC with a reference channel 
and post-TPD estimation 
& subtraction 


strong 
RFI 


impulse-like 


Removing (blanking) 
in temporal domain 


Removing (blanking) 
in temporal domain 


narrow- band 


Removing (blanking) 
in spectral domain 


ANC with a reference 
channel; HOS analysis 


wideband 


ANC with 

a reference channel 


ANC with 

a reference channel 


Connected 
radio 

interferometer 


weak 
RFI 


impulse- like 


tolerable 


tolerable 


narrow- band 


tolerable 


Post-correlation ANC with 
reference channels; HOS analysis 


wideband 


tolerable 


Post correlation ANC with 
reference channels 


strong 
RFI 


impulse- like 


Removing (blanking) 
in temporal domain 


Removing (blanking) 
in temporal domain 


narrow-band 


Removing (blanking) 
in spectral domain 


ANC with a reference 
channel. Spatial filtering 


wideband 


ANC with a reference 
channel. Spatial filtering 


ANC with a reference 
channel. Spatial filtering 


VLBI 


weak 
RFI 


impulse-like 


tolerable 


tolerable 


narrow-band 


tolerable 


tolerable 


wideband 


tolerable 


tolerable 


strong 
RFI 


impulse- like 


Removing (blanking) 
in temporal domain 


Removing (blanking) 
in temporal domain 


narrow-band 


Removing (blanking) 
in spectral domain 


ANC with a reference 
channel 


wideband 


ANC with a reference 
channel 


ANC with a reference 
channel 



RFI - radio frequency interference. 
TPD - total power detector. 
ANC - adaptive noise cancellation. 
HOS - higher order statistics. 
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List of captions 

Fig. 1. a) The empirical influence function as a function of RFI amplitude, cr = 1, n = 10^, the estimates 
are made with initial data. The curves here and elsewhere are parameterized by the fraction of RFI in the 
total volume of data: e = 0.01; 0.025; 0.5; 

b) estimates made using trimmed data; 

c) estimates made using winsorized data; 

d) estimates made as the median of pairwise averaged squares', 

e) estimates made using the Qn estimate algorithm; 

f) estimates made using the biweight variance algorithm; 

g) estimates made using the bend midvariance algorithm; 

h) estimates made using the algorithm of exponential weighting-, 

i) estimates made using the algorithm of interquartile range. 

Fig. 2. Block diagram of computer simulations: winsorization is used for RFI mitigation when the total 
power detector is the backend output as in section 4. 1. 

Fig. 3. Results of computer simulations with the algorithm shown in Fig. 2: 

a) noise with the normal distribution, /i = 0.0, cr = 0.5, no interference; 

b) total power detector output, each point in this figure corresponds to squaring and averaging of 10^ samples 
in figure a), there are two steps, "up" at point #100 and "down" at point #200 corresponding to the increase 
of a from value 0.5 to the value a + Act, Act = 0.05. 

c) interference is added to the noise a): random impulses with the Poisson distibution (A = 0.04) and the 
lognormal disrtibution of amplitudes (mean=10, standard deviation=5); 

d) total power detector output with input signal c); 

e) total power detector output with input signal c) and preliminary winsorization (equation (9)), note the 
difference of scale in d) and e). 

Fig. 4. Block diagram of computer simulations: exponential weighting is used for RFI mitigation in the 
frequency domain when the total power detector is the backend output as in section 4. 2 (pulsar observations) 
or the correlator as in section 4. 3. (radiointerferometric observations). 

Fig. 5. Results of computer simulations of exponential weighting with the algorithm shown in Fig. 4: 

a) time- frequency presentation of power spectrum consisting of system noise, emission and absorption lines 
and RFI (randomly binary-phase manipulated signals), L=50 time sections of spectrum divided on 256 
channels, each spectrum is the mean of M = 100 instantaneous spectra; 

b) power spectrum obtained as a solution to equation (18) for each spectral channel and M = 100 samples; 

c) power spectrum averaged in time, using L = 50 sample spectra from a); 

d) power spectrum averaged in time, using L = 50 sample spectra from b), i.e., with RFI mitigation. 
Fig. 6. Results of computer simulations of exponential weighting with the algorithm shown in Fig. 4: 

spectra of interference overlap with spectral lines. Other parameters are similar to Fig. 5, 
a) time- frequency presentation of the power spectrum consisting of system noise, emission and absorption 
lines and RFI; b) power spectrum obtained as a solution to equation (18) for each spectral channel and 
M 100 samples; 

c) power spectrum averaged in time, using L = 50 sample spectra from a); 

d) power spectrum averaged in time, using L = 50 sample spectra from b), i.e., with RFI mitigation. Both 
spectral lines are clearly visible. 

Fig. 7. Results of computer simulations of exponential weighting with the algorithm shown in Fig. 4 
when "dirty" and "clean" signals are applied to the correlator; L=50 time sections of the spectrum divided 
on 256 spectral channels, each spectrum is the mean of M = 100 instantaneous spectra; the spectra at the 
first input of the correlator are shown, the spectra at the second input are similar to a) and b): 

a) time- frequency presentation of power spectrum consisting of system noise, and RFI (frequency- modulated 
bursts) ; 

b) power spectrum after RFI mitigation - exponential weighting of each instantaneous spectrum using vari- 
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ances obtained as a solution to equation (18) for each spectral channel and M = 100 samples; 

c) cross-correlation in the presence of RFI and no RFI mitigation, the central 200 channels are shown; 

d) cross-correlation in the presence of RFI and with RFI mitigation; 
take notice of the difference of the vertical scales in c) and d). 

Fig. 8. Examples of RFI mitigation at the Effelsberg radio telescope during observations in continuum 
at central frequency 1645 MHz, bandwidth 20MHz. A selection of eight scans of the source 1448+762 is 
represented. Pairwise records were made simultaneously for the channel with RFI mitigation and the channel 
without RFI mitigation. 

Fig. 9. Radio image of the source 1448+762 built using scans similar to those in Fig. 2: left panel - 
without RFI mitigation, right panel - with RFI mitigation. 

Fig. 10. Pulsar B0329+54.07 observed at WSRT at 1625MHz. Data were recorded during 10 sec, 40 
Msamples/sec. Upper row, left panel: TPD ouputs for two polarizations, raw data with RFI, right panel: 
TPD output, RFI removed. Middle row, left panel: time fragment of the running power spectrum with RFI; 
the same time fragment, RFI removed. Lower row. Left panel: pulsar profile made with raw data over 10 
sec, middle panel: pulsar profile, RFI removed, right panel: pulsar profile observed at 1420 MHz, no RFI. 

Fig. 11. Radio images of the source DA240 observed at WSRT, central frequency 357 MHz, bandwidth 
20 MHz. Upper row, left panel: image without RFI mitigation; right panel: image with RFI mitigation. 
Lower row: central parts of the images with and without RFI. 

Fig. 12 . Radio images of the source 4C34.47 observed at WSRT, bandwidth 20 MHz. Upper row, left 
panel: central frequency 345 MHz, without RFI mitigation system; right panel: central frequency 345 MHz, 
with RFI mitigation system. Notice the difference of the intensity levels in the figures. Low row ("toxicity" 
test), left panel: central frequency 1420 MHz without RFI mitigation system, right panel: central frequency 
1420 MHz with RFI mitigation system. After subtracting one image from the other the rms noise is less 
than 0.7 mJy/beam which signifies a good similarity in the images. 
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Fig. 1. — a) The empirical influence function as a function of RFI amplitude, a = l,n = 10^, estimates 
made with initial data. The curves here and elsewhere are parameterized by the fraction of RFI in the 
total volume of data:e = 0.01; 0.025; 0.5; b) estimates made using trimmed data; c) estimates made using 
winsorized data; d) estimates made as the median of pairwise averaged squares ; e) estimates made using the 
Qn estimate algorithm; f) estimates made using the biweight variance algorithm; g) estimates made using 
the bend midvariance algorithm; h) estimates made using the algorithm of exponential weighting . 
i) estimates made using the algorithm of interquartile range. 
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Fig. 2. — Block diagram of computer simulations: winsorization is used for RFI mitigation when the total 
power detector is the backend output as in section 4. 1. 
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Fig. 3. — Results of computer simulations with the algorithm shown in Fig. 2: a) noise with the normal 
distribution, = 0.,(j = 0.5, no interference; b) total power detector output, each point in this figure 
corresponds to squaring and averaging of 10^ samples in figure a), there are two steps, "up" at point #100 
and "down" at point #200 corresponding to the increase of a from value 0.5 to the value (7 + Act, Act = 0.05. 
c) interference is added to the noise a): random impulses with the Poisson distibution (A = 0.004) and the 
lognormal disrtibution of amplitudes (mean=12, standard deviation=6); d) total power detector output with 
input signal c); e) total power detector output with input signal c) and preliminary winsorization, note the 
difference of scale in d) and e) . 
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Fig. 4. — Block diagram of computer simulations: exponential weighting is used for RFI mitigation in the 
frequency domain when the total power detector is the backend output as in section 4. 2 (pulsar observations) 
or the correlator as in section 4. 3. (radiointerferometric observations). 
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Fig. 5. — Results of computer simulations of exponential weighting with the algorithm shown in Fig. 4: a) 
time-frequency presentation of power spectrum consisting of system noise, emission and absorption lines and 
RFI (randomly binary-phase manipulated signals), L=50 time sections of spectrum divided on 256 channels, 
each spectrum is the mean of M = 100 instantaneous spectra; b) power spectrum obtained as a solution to 
equation (18) for each spectral channel and M = IcRPsamples; c) power spectrum averaged in time, using 
L = 50 sample spectra from a); d) power spectrum averaged in time, using L = 50 sample spectra from b), 
i.e., with RFI mitigation. 




Fig. 6. — Results of computer simulations of exponential weighting with the algorithm shown in Fig. 4: 
spectra of interference overlap with spectral lines. Other parameters are similar to Fig. 5, a) time- frequency 
presentation of the power spectrum consisting of system noise, emission and absorption lines and RFI; b) 
power spectrum obtained as a solution to equation (18) for each spectral channel and M = 100 samples; 
c) power spectrum averaged in time, using L = 50 sSiiple spectra from a); d) power spectrum averaged in 
time, using L = 50 sample spectra from b), i.e., with RFI mitigation. Both spectral lines are clearly visible. 




Fig. 7. — Results of computer simulations of exponential weighting with the algorithm shown in Fig. 4 when 
"dirty" and "clean" signals are applied to the correlator; L=50 time sections of the spectrum divided on 
256 spectral channels, each spectrum is the mean of M = 100 instantaneous spectra; the spectra at the first 
input of the correlator are shown, the spectra at the second input are similar to a) and b): a) time- frequency 
presentation of the power spectrum consisting of system noise, and RFI (frequency- modulated bursts); b) 
power spectrum after RFI mitigation - exponential weighting of each instantaneous spectrum using variances 
obtained as a solution to equation (18) for each spectral channel and M = 100 samples; c) cross-correlation 
in the presence of RFI and no RFI mitigation, the central 200 channels are shown; d) cross-correlation in 
the presence of RFI and with RFI mitigation; take notice of the difference of the vertical scales in c) and d). 
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Fig. 8. — Examples of RFI mitigation at the Effelsberg radio telescope during observations in continuum 
at central frequency 1645 MHz, bandwidth 20MHz. A selection of eight scans of the source 1448+762 is 
represented. The pairwise records were made simultaneously for the channel with RFI mitigation and the 
channel without RFI mitigation. 
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Fig. 9. — Radio image of the source 1448+762 built using scans similar to those in Fig. 2: left panel 
without RFI mitigation, right panel - with RFI mitigation. 



24 



Pul ^r B0329+54.07, 1625MHz, 29June 2005, average int en/a 1=0. 001 sec 




Pulsar B0329-f54.07, 1625MHz, 29June 2005, average interval=0.001 sec 



1 . , 1 1 









□ 1 23456789 10 

sec 





1 1 

.. .J \ 


1 1 

. j 1 .. 


1 1 1 1 

1 1 1 1 


,.J. ... 





23456789 10 
sec 



Pulsar B0329454.07, 1625MHz, 29June 2005, 
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Fig. 10.— Pulsar B0329+54.07 observed at WSRT at 1625MHz. Data were recorded during 10 sec, 40 
Msamples/sec. Upper row, left panel: TPD ouputs for two polarizations, raw data with RFI, right panel: 
TPD outputs, RFI removed. Middle row, right panel: time fragment of the running power spectrum with 
RFI; the same time fragment, RFI removed. Lower row, left panel: pulsar profile made with raw data over 
10 sec, middle panel: pulsar profile, RFI removed, right panel: pulsar observed at 1420 MHz, no RFI. 
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Fig. 11. — Radio images of the source DA240 observed at WSRT, central frequency 357 MHz, bandwidth 20 
MHz, Upper row, left panel: image without RFI mitigation; right panel: image with RFI mitigation. Lower 
row: central parts of the images with and without RFI. 
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Fig. 12. — Radio images of the source 4C34.47 observed at WSRT, bandwidth 20 MHz. Upper row, left 
panel: central frequency 345 MHz, without RFI mitigation system; right panel: central frequency 345 MHz, 
with RFI mitigation system. Notice the difference of intensity levels in the figures. Lower row, left panel: 
central frequency 1420 MHz without RFI mitigation system, right panel: central frequency 1420 MHz with 
RFI mitigation system. After subtracting one image from the other the rms noise is less than 0.7 mJy/beam 
which signifies a good similarity in the images. 
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